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Least squares (LS) fitting is one of the most fundamental techniques in science and engineering. It 
is used to estimate parameters from multiple noisy observations. In many problems the parameters are 
known a-priori to be bounded integer valued, or they come from a finite set of values on an arbitrary 
finite lattice. In this case finding the closest vector becomes NP-Hard problem. In this paper we propose a 
novel algorithm, the Tomographic Least Squares Decoder (TLSD), that not only solves the ILS problem, 
better than other sub-optimal techniques, but also is capable of providing the a-posteriori probability 
distribution for each element in the solution vector. The algorithm is based on reconstruction of the 
vector from multiple two-dimensional projections. The projections are carefully chosen to provide low 
, computational complexity. Unlike other iterative techniques, such as the belief propagation, the proposed 

algorithm has ensured convergence. We also provide simulated experiments comparing the algorithm to 



(N 



O 

o 



other sub-optimal algorithms. 



Index Terms 



Integer Least Squares, Bayesian decoding, sparse linear equations. MIMO communication systems. 

I. Introduction 

A multiple-input-multiple-output (MIMO) system is a communication system with d transmit antennas 
and p receive antennas. The tap gain from transmit antenna i to receive antenna j is denoted by Hy. In 
each use of the MIMO channel a signal vector s = (s%, Sd) r is independently selected from a set of 
constellation points A according to the data to be transmitted, so that s £ A d . The received vector x is 
given by: 

x = Hs + n (1) 
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The vector n is an additive noise in which the noise components are assumed as zero mean, statistically 
independent Gaussians with a known variance a 2 . The channel matrix which is assumed to be known, 
comprises i.i.d. elements drawn from a (circularly symmetric zero-mean complex) normal distribution of 
unit variance. In the case where the MIMO linear system is complex-valued we use the standard method 
to translate it into an equivalent double-size real-valued representation that is obtained by considering the 
real and imaginary parts separately. The MIMO detection problem is then becomes finding the transmitted 
vector s given H and x. The optimal maximum likelihood (ML) solution is: 

s = arg min ||Hs — xll 2 (2) 

However, ML decoding has exponential computational complexity which makes it unfeasible when either 
the number of transmitted antennas or the constellation size are large. Actually, for a general H, it is 
known to be NP-HARD both in the worst-case sense (H as well as in the average sense |2|. It can be 
easily verified that the MIMO ML detection problem is equivalent to a least square lattice search problem 
that is known to be NP hard. A simple approximation is the zero-forcing (ZF) algorithm which is based 
on a linear decision ignoring the finite constellation constraint: 

s = (B T BQ -1 H T x (3) 

and then, neglecting the correlation between the symbols, finding the closest constellation point for each 
symbol independently. This scheme performs poorly due to its inability to handle ill-conditioned channel 
matrix realizations. Somewhat better performance can be obtained by using a minimum mean square 
error (MMSE) filter instead of ZF on the un-constrained linear system: 

s= (H T H + o- 2 /) _1 H T x (4) 

and then finding the closest lattice point in each component independently. A vast improvement over 
the linear approach can be achieved by using sequential decoding. This algorithm, known as MMSE 
V-BLAST or MMSE-SIC, has the best performance for this family of linear-based algorithms. However, 
there is a still a significant gap between the detection performance of the MMSE-SIC algorithm and the 
performance of the optimal ML detector. The complexity of all these algorithms is 0(p 3 ) where p is the 
number of receive antennas (we assume p > d). These algorithms can also easily provide probabilistic 
(soft-decision) estimates for each symbol (or each bit). 

Many alternative structures have been proposed to approach the ML detection performance. For 
example, the sphere decoding algorithm (3]> approaches using the sequential Monte Carlo framework 
H and methods based on semidefinite relaxation |5), (6l have successfully been implemented. Although 
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the detection schemes listed above have significantly reduced computational complexity, sphere decoding 
is still exponential in the average case Q) and semidefinite relaxation is high-degree polynomial. Neither 
of these approaches can be easily used in real-world hardware architecture applications. Since these algo- 
rithms find the closest point in the lattice, it is not straight-forward to compute a-posteriori probabilities 
per symbol or per bit, which is required, when forward error correction is used (e.g., in communication 
applications), but it can be done with increased complexity of the sphere decoding procedure [8 ]. Thus, 
there is still a need for low complexity detection algorithms that can achieve good performance with 
low-order polynomial complexity, that are capable of providing per-bit likelihood ratios. 

In this paper we propose a novel iterative technique, which we dub Tomographic Least Squares Decoder 
(TLSD), that is based on 2D projections followed by iterative optimization. The solution also allows us 
to provide a-posteriori probability distributions to each bit of each variable, something desirable in coded 
communication systems. Such probabilities are more complicated to evaluate using sphere decoding types 
of solution. 

The paper proceeds as follows. In Section 2 we present the proposed TLSD algorithm. Experimental 
results are shown and discussed in Section 3. 

II. Tomographic Decoding of Constrained Linear Systems 

In this section we present a novel polynomial time algorithm for solving the bounded integer least 
squares problem. The algorithm outperforms other reduced complexity algorithms with lower complexity. 
The algorithm has two important steps. The first step is translating the problem into a set of two- 
dimensional problems. The second step comprises of solving iteratively the two-dimensional problems 
by using data received from other two-dimensional problems. This is very similar to tomographic imaging, 
where an object is reconstructed from its projections on lower-dimensional subspaces. Hence we dub it 
Tomographic Least Squares Decoder (TLSD). The difference is that our object is discrete, and the data 
that is shared among the projections consists of probability distributions. The second step can also be 
interpreted as an instance of the incremental EM algorithm J9). This will allow us to prove the convergence 
of the algorithm. 

A. The two-dimensional projections 

Our approach can be viewed as a combination of a two-dimensional generalization of the ZF solution 
with optimal solution of the 2D subsystems obtained by this generalization. 
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Let hi, hd be the columns of H and for each 1 <i< j <d let Ay be the matrix obtained from H 
by removing both i-th and j-th columns. It can be easily verified that the transformation: 

Py = I — Aij (A^ Aij) Ay (5) 

is an orthogonal projection into the complement of the sub-space spanned by {hk\k ^ i,j}- Hence 

PyHs = Py ^2 hkSk = '''V 1 ''-'' + P 'Ar s V ( 6 ) 

k 

Applying the linear transformation Py on both sides of the equation Hs + n = x, yields a set of p 
equations that depends only on the two variables Sj and Sj\ 

PijhiSi + PijhjSj + Pyn = Pyx (7) 

Using the simplifying notation: Hy = Py[hj,hj], ny = Pyn and xy = Pyx, Eq. © can be written as: 

Hy[sj, Sj] T + ny = Xy (8) 



where ny ~ jV(0, cx 2 Py). The density function of 



fij(^ij] Si, Sj) — 2 6Xp( 2 || x j.j Hy[sj,Sj] || ) (9) 



Note that this is a two dimensional density function since the vector Xy belongs to a two dimensional 
subspace spanned by Pyh, and Pyhj. Furthermore, the orthogonal projection of an isotropic Gaussian 
variable is still isotropic in the projected space. 

We have converted the original linear system into ( 2 ) systems of sparse linear equations. If we take 
only non-overlapping projections (e.g. Pi2,P34, •••,Pd-id) and solve the corresponding linear systems, it 
can be easily verified that we get exactly the linear ZF solution. Our approach is based on taking all the 
(2) sparse systems. Due to the overlap between the projections, each of the solvers of the sub-problems 
provides information to the other solvers. 

Ignoring the noise correlation between equation sets obtained by different projections, the likelihood 
function of s € A d , based on the sparse linear systems: 

Hij[si,Sj] + ny = Xy , l<i<j<d (10) 



is: 



/(x;s) = Y[fij(^ij)Si,Sj) (11) 

i<j 

^ 6XP( " 2^ XJ H x ^' " H 'J k' I' 2 ) 



5 



Note that /(x;s) is not the precise likelihood function since we ignore the noise correlation between 
equations derived from different projections. Note however that all pairwise correlations are still captured 
by the relevant 2-D subproblems, so basically we only give up noise correlations of order 3 and higher. 

Our goal now is finding the maximum-likelihood solution of the new system: s = argmax s /(x;s). 
In the next section we present an iterative method for maximizing /(x; s). The main point of this paper 
is that by applying the 2D projections we shift from the original likelihood function into a very similar 
function that is much easier to optimize. The sparsity of the new system makes /(x; s) a much smoother 
function than the original likelihood function. This smoothness enables applying an effective iterative 
search. A similar situation occurs in LDPC codes[10] where the sparsity of the parity-check matrix 
results in a smooth likelihood function. 

B. Iterative solution of the sparse problem 

We have now converted the original linear system, into (fj sets of sparse equations. The second step 
of our approach comprises of solving iteratively the two-dimensional problems by using data received 
from other two-dimensional problems. 

Given an a-priori probability vector for Sj , Sj we can now easily use Xy to update these probabilities 
in a locally optimal way. Assume that for each i = 1, d we have an a-priori probability distribution on 
Si i.e., probability vectors 6i = 0j(M)), where 9i{k) = p (sj = a k ), where A = {ai, om} is 

the finite symbol set. Given x^ we can compute the a-posteriori probability for Sj , sj denoted by 0\ , 6°j 
respectively and given by 

M 

6f(k) oc Oiik^O^Dijia^ai) (12) 

M 
k=l 

where 

Dij(a k , at) = fij(xij;a k , at) (13) 

and the notation oc indicates normalization of the vector to make it a distribution. We can now iterate 
the updates of Oi,i = l,...,d, by choosing at each iteration a new pair i < j and updating 6i,6j. It can 
be shown that this is an instance of the EM algorithm. 
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To initialize the process we need a good choice of the a-priori probability vectors 0j. This can be 
done for example using a soft version of the ZF solution. Let P, be the ZF one-dimensional orthogonal 
projection into complement of the subspace spanned by {h^jA; ^ i}. Then the initial parameter values 
are: 

Oi(k) oc exp(-^ HPihiOfc - Pixf) (14) 

C. The TLSD algorithm 

To decode an integer LS problem we perform the following: We first compute all the matrices H^. 
This amounts to (f\ QR factorizations for each i < j. This has complexity O (<i 2 p 3 ), but it is done once 
in the beginning of the decoding process. Now for each received vector x we first use a ZF receiver to 
generate the a-priori probability distributions 6i,i = 1, . . . ,d. This has complexity of 0(p 3 ), since the 
main problem is the computation of the ZF receiver. Computing the priors is 0{Md). 

The next step is to compute for each two-dimensional vector of constellation points the metric using 
equation ( TT3T ). This has complexity of O (<i 2 M 2 ). This is done once for each received vector. Now we go 
over all vectors Xjj sequentially and update 9i,6j using (fl2l ). This is done until convergence is achieved, 
typically with few iterations. The overall complexity is 0(M 2 Nu eT ). After convergence the we obtain 
a-posteriori probabilities per symbol. A hard-decision solution is given by choosing for each i = l,...,d 
the most probable symbol: 

§i = arg max OAk) (15) 
l<fc<M 

The algorithm-box in Table I summarizes the TLSD algorithm. The proposed TLSD algorithm is based 
only on two-dimensional subspaces. It is straightforward to improve the algorithm by using projections on 
higher dimensional subspace. This can have improved performance, and higher computational complexity. 
Finally we compare the likelihood of the solution vector with the likelihood of the MMSE-SIC solution 
and choose the one with higher likelihood. It turns out that since these algorithms use different type of 
information about the solution, that this improves the performance, especially for low SNR situations. 

III. Simulations 

In this section we provide simulation results for the proposed detector over various uncoded MIMO 
systems. We assume a quasi-static fading channel with a frame length of 100. Under the assumption 
of block-fading channel model, the channel matrix H is constant for 100 channel uses. The channel 
matrix comprised i.i.d. elements drawn from a zero-mean normal distribution of unit variance. We have 
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TABLE I 

The Tomographic Least Squares Decoder (TLSD) 



Input: An integer LS problem: Hs + n = x, a noise level a 2 
and a finite symbol set {ai, »m}. 

Initialization: 
For i = 1, d 

initialize (8i(k) : a k £ A) using zero-forcing d!4t . 

For each pair 1 < i < j < d 

Compute the projection Pjj using Eq. ((5} and compute: 

Dij(a k ,ai) = exp(- 2^2-||P 4j (x - hid k - hja e )\\ 2 ) ■ 

End 



Do until convergence 
For 1 < i < j < d 

Update the distributions 9i,9j\ 

Oi{k) <x 0i(k) E"i ^(OAi(o fc ,ai) 

End 
End 

For i = 1, d 

Choose Si = argmaxfc 9i(k). 



used 10,000 realizations of channel matrix. This results in 10 vector messages. The SNR is defined as 
10 log 10 (£'b/A r o) where Ef, is the average received energy per symbol at each receiver antenna. 

Fig. 1 shows the symbol error rate (SER) versus SNR for a 8 x 8 BPSK MIMO system. The performance 
of the TLSD method is compared to ML detection and to other linear suboptimal algorithms: the linear 
MMSE and the sequential MMSE V-BLAST. In all our experiments the number of TLSD iterations was 
limited to 10. It can be seen that the TLSD algorithm is significantly better than the MMSE-SIC at 
the same computational complexity. Fig. 2 depicts similar results for a 16 x 16 4-PAM MIMO system. 
The TLSD decoder significantly outperforms the MMSE-SIC, while having comparable computational 
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complexity. 

IV. Conclusions 

Solving integer least squares problems is an important problem in many fields. We have proposed a 
novel technique based on tomographic principle of reconstruction from projections. We showed that the 
method always converges. Furthermore, the proposed method has good performance competitive to all 
other polynomial algorithms for solving the problem as demonstrated in simulations. Finally the method 
can be extended to provide a-posteriori probabilities per bit for use in coded communication systems or 
combined with sphere decoding, to improve its performance. 
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Fig. 1. Results for 8 x 8 system, A = {-1, 1} 
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